Electronic properties of graphene multilayers 
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We study the effects of disorder in the electronic properties of graphene multilayers, with special 
focus on the bilayer and the infinite stack. At low energies and long wavelengths, the electronic self- 
energies and density of states exhibit behavior with divergences near half-filling. As a consequence, 
the spectral functions and conductivities do not follow Landau's Fermi liquid theory. In particular, 
we show that the quasiparticle decay rate has a minimum as a function of energy, there is a universal 
minimum value for the in-plane conductivity of order e^ /h per plane and, unexpectedly, the c-axis 
conductivity is enhanced by disorder at low doping, leading to an enormous conductivity anisotropy 
at low temperatures. 

PACS numbers: 81.05.Uw 73. 21. Ac 71.23.-k 



Introduction. Recently, many properties of a single 
graphene sheet have been studied theoretically by sev- 
eral groups. These properties are, in many cases, found 
to be unconventional due to the peculiar band structure 
of graphene which is described in terms of Dirac fermions 
at the edge of the Brillouin zone (BZ). This activity was 
triggered by the realization that graphene could be ob- 
tained and studied experimentally 1], and the subse- 
quent experiments that followed to further characterize 
the material |3,l^- 

More recently the attention has turned to graphene 
multilayers 4] and, in particular, to bilayers that also 
show anomalies in the integer quantum Hall effect 
(IQHE) 5] and have received theoretical attention ja^Q- 
In this paper we show that the bilayer graphene has also 
very unconventional behavior in its spectral and trans- 
port properties. These properties, however, are quite dif- 
ferent from those of the single layer. In fact, the anoma- 
lous behavior is a property of all multilayer graphene 
systems we have studied. It is also worth noting that 
the higher complexity of the multilayer systems may be 
of interest for device applications since it also allows for 
greater flexibility in tailoring the electronic properties. 

There are two key ingredients that make the physics 
of the graphene multilayers unconventional. Firstly, due 
to the relatively weak interlayer coupling it inherits some 
properties from its parent material, the single graphene 
sheet. The existence of Dirac points in the BZ, where 
the electron and hole bands becomes degenerate, arises 
from the physics of the single layer in conjunction with 
the second important ingredient that is the peculiar ge- 
ometry that results from the A-B stacking of the planes. 
This geometry implies that the binding of the planes due 
to orbital overlap is mainly sitting on one of the sublat- 
tices (that we call A, the other sublattice being B) in 
each plane (the different planes of the units cell are de- 
noted by 1 and 2). The main residence of the low-energy 



states is then on the sublattice B. Nevertheless, electron 
transport coming from nearest neighbor hopping must 
go over the A sublattice. This feature implies that even 
though the total density of states at half-filling is finite, 
the density of states on the A sublattice vanishes as the 
Dirac point is approached, leading to unconventional in- 
plane and out-of-plane transport properties. The main 
purpose of this work is to show how this comes about 
and to highlight the unusual behavior of the self-energies 
due to disorder near half-filling. 

Graphene bilayer. At low energies and long wave- 
lengths, one can expand the electronic spectrum close 
to the K and K' points in the BZ, leading to a Hamilto- 
nian of the form % H = J^j^ ^t(k)Ho(k)^(k), where 
k = {kx^ ky) is the two-dimensional momentum measured 
relative to the K {K') point. 
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^^k) = (c^, k 4i,k CA2,k 42,k) is the electron spinor 
creation operator, (/)(k) = idin~^ (ky / kx) is the two- 
dimensional angle in momentum space, vy = 3ta/2 is 
the Fermi-Dirac velocity, t ^ 3 eV is the in-plane hopping 
energy, a ~ 1.42 A is the interatomic distance within the 
layers, t± ~ 0.35 eV is the interlayer hopping energy p. 
The 2 X 2-blocks on the diagonal of (^ are identical to 
the continuous approximation that leads to the massless 
Dirac spectrum in a single layer. In what follows we use 
units such that v^ = I = h, and suppress the spin and 
valley indices unless otherwise specified. The Hamilto- 
nian Q can be easily diagonalized and the energy spec- 
trum is made out of four bands with energy: t±/2^E{k) 
and -t±/2 ± E{k), where E{k) = ^/t\jtn^. The re- 
sulting spectrum is made out of two vertex touching hy- 
perbolae at zero energy, separated by a gap of energy t±_ 



from two other hyperbolae. 

In what fohows it is convenient to introduce a local 
frequency dependent self-energy which, due to the lattice 
structure, has different values in the A (11^(0;)) and the 
B {T^b{^)) sublattices. We take these into account by a 
diagonal matrix H^{uj) and the Green's function is then 
given by: 



G-^(cj,k) 



Ho{k)-H^{uj). 



Because the Hamiltonian factorizes into two blocks it is 
simple to work out the explicit expression for the Green's 
function. G(a;, k) is a 4 x 4 matrix, but for our purposes 
here the most important components are the diagonal 
ones that are given by: 



Effects of disorder. We use standard techniques to 
study the effect of disorder and average over impurity 
positions to get the disorder-averaged propagators flOJ . 
The effect of repeated scattering from a single impurity 
can be encoded in a self-energy which can be computed 
in the T- matrix approximation as: 



S(cj) = V[l-G{u;)V]-^/N, 
GH = ^G(c^,k)/7V, 



(3) 



where V is the strength of the impurity potential, N is 
the number of units cells in each plane. We turn the 
momentum sum into an integral by the introduction of 
a cut-off, A (^ l/<^), that we estimate (A ^ 7eV) by a 
Debye approximation conserving the number of states in 
the BZ 7] . Due to the special form of the propagator and 
the impurity potential, the self-energy is diagonal. The 
result for unitary scattering (or site dilution) is obtained 
by introducing a finite density of vacancies ni and taking 
the limit F ^ oo, and one finds: 

Sa(cj) = -ni/GAA(^) , Sb(^) = -ni/GBB(^)- (4) 

In the full Born approximation (FBA) one uses the bare 
propagators on the right hand side of Q), the resulting 
self-energies are linear in the impurity concentration so 
that this approximation neglects all correlations between 
different scattering centers. In the coherent potential ap- 
proximation (CPA) |jj| one assumes that the electrons 
are moving in an effective medium that is characterized 
by the self-energies Ha and Hb- These must be deter- 
mined self-consistently by using the dressed propagators 
on the right hand side of Q), thus this approximation 



includes some effects of correlations between the scatter- 
ing events (it does not describe Anderson localization). 
Using the explicit form of the propagators we obtain the 
following self-consistent equations: 
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(5) 



These equations include intervalley scattering in the in- 
termediate states. It is also straightforward to obtain the 
density of states (DOS) on sublattice a {a = A^ B) from 

Pa(^) = — Ii^G^^/tt. In the clean case, one gets: 
pliio) = \io\[l + Q{\io\-t^)]/2A\ 

where Q{x) is the Heaviside step function. Notice that 
the DOS on the A sublattice vanishes at zero energy (as 
in the case of a single layer), while Pb(0) is finite. 

We have solved the CPA equations in Eq. Q and 
an example of the resulting self-energies and the corre- 
sponding DOS are shown in Fig. Q Within the FBA 
E^ diverges as n\K^ /uo up to logarithmic corrections as 
a; ^^ 0. In the single layer, the CPA makes the self-energy 
finite at a; = 0. In contrast, the bilayer (and the multi- 
layer) does not allow a finite E^ at cj = in the CPA. 
One can see this by studying Eq. Q at cj = 0. Then the 
last line of Eq. ((SJ divided by E^ must vanish, and this is 
not possible for finite values of E^, furthermore E^ must 
vanish in this limit. These results imply that, contrary 
to the single plane case, pa(^ -^ 0) is zero even within 
the CPA. More explicitly, by defining E^iE^ = — ^A^ one 
can show (assuming E^ ;» t±_ and E^ ^ oo) that T^a and 
Eb are given asymptotically by: 



Ea(c^) 



r/3 {tl^^K^ln, 
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-1/3 



Eb(c.) ^ e-^'-/\n-Xi/tlf'uj'^\ 



(6a) 
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for u <C A^^^/nit^, where ^ is given by n\ = ^ln(l/^). 
Hence, the self-energies are rather anomalous and clearly 
violate Landau's Fermi liquid theory. Note that V^A ~ 
y^A is exactly the scale generated by disorder in the 
single plane case [l^ . 

Infinite stack. The extension of the bilayer model to 
multilayers is straightforward. Upon a Fourier transfor- 
mation in the c-direction we can immediately use the the 
Hamiltonian in Eq. Q with k = {kx^ky^k±)^ where k± 
is the momentum along the c-axis (— 7r/2 < k±d < 7r/2), 
and make the substitution t± -^ 2t±cos{k±d)^ where 
d ^ 2.5a is the interlayer distance |3|. The calculations 
of the DOS and self-energies proceed as previously but 
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FIG. 1: (Color online) Bilayer DOS and self-energy (in units 
of the A) in the CPA approximation. DOS on the A (a) and 
B(c) sublattice as a function of the frequency (in units of 
A), imaginary part of the self-energy on the A(b) and B(d) 
sublattice as a function of the frequency. 




0.05 0.1 0.15 0.05 0.1 0.15 0.05 0.1 0.15 

vpkw/A vpkii/A vphii/A 

FIG. 2: Intensity plot of the multilayer spectral function (0 
in the u (in units of A) versus in-plane momentum k\\ (in 
units of A/i;F) plane, for different values of the transverse 
momentum k± (in units of 1/d) and impurity concentration. 
Upper row: m — lO"'^; lower row: m = 5 x 10~^. 



with one additional /c^-integral 12]. The electron spec- 
tral function, that is measurable in angle resolved photo- 
emission (ARPES), is given by: 

A(k,u;) = -(2/7r)Im[G?A(k,t^) +GgB(k,a;)]. (7) 

In Fig. Owe present an intensity plot of the spectral func- 
tion for three values of the perpendicular momentum and 
two impurity concentrations. In Fig.lHJwe show four con- 
stant k-cuts. It is clear that disorder leads to broad peaks 
and that the high-energy branch is less affected by the 
disorder than the low-energy branch. Electron-electron 
interactions lead to an extra contribution to the self- 
energy (not included in the plots) that scales linearly with 
frequency 13], ImSeel^) oc |a;|. Hence, disorder leads to 
a quasiparticle damping that increases at low frequencies, 
while the electron-electron contribution increases at high 
frequencies, the final result is a quasiparticle decay rate 
that has a minimum at a finite frequency. This highly 
non-Fermi liquid behavior was also found in the case of 
single layer graphene 10] and is present in all multilayer 
systems. Notice also that the peak positions are shifted 
by disorder. Another interesting feature is the appear- 
ance of a new peak in the spectrum near k±d = 7r/2, 
which is clearly visible in Fig. IHlJc) for ni = 10~^. This 
extra peak is due to the fact that the real part of the 
self-energies can act as a "mass" term in the Dirac equa- 
tion leading to the formation of a "pseudogap" . For even 
higher values of the disorder the peak goes away and only 
a broad shoulder remains. Once again, we can clearly see 
violations of Fermi liquid theory. 

Transport properties. The conductivity can be com- 
puted from the Kubo formula: 

a{u;) =1 dte'^\[J{t), J{^)]) /{Suj)=i]l{uj)/{uj + i5), (8) 
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FIG. 3: (Color online) Constant k cuts of the spectral func- 
tion in the multilayer as a function of the frequency (in units 
of A), (a) k± = 0, /cii = 0, (b) k± = 0, /c,, = 0.05A, (c) 

k± = 7r/2, k\\ = 0, (d) k± = 7r/2, k\\ = 0.05A. 



where J is appropriate the current operator, and 11 (cj) 
the associated current-current correlation function. The 
current operators are dictated by gauge invariance [l^ 
and are computed using the Peierls substitution. For 
instance, for plane 1 with the current in the x di- 
rection we find: Jxi = -^^pe X]k[^*''^^<^Ai,k^Bi,k ~ 
g-z7r/3^-^ ^^^^ ^j . £q^ |-]^g multilayer, the c-axis cur- 
rent operator is: J± = —2et±d'^^sm{k±)[c]^ j^c^ j^ + 

^A k^Ai k] • ^^ evaluate the current-current correlators 
with the disorder- averages propagators in the Matsub- 
ara formalism 15]. For the real part of the frequency- 
dependent conductivity we find: 



Re [cr{uj)] /ao= de{nY{uj) - nF(cc; + e))S(e, e - 



-u;)/uj, (9) 



where nF(e) is the Fermi-Dirac function, ctq is the unit 
of conductance (cro|| = 4e^/7r/i for the bilayer, and cro|| = 
c'"o±('^F/2tx<i)^ for the multilayer), S is a kernel that we 
have evaluated [l^ . Note that since we model the impu- 
rities as purely local, there are no vertex corrections. We 
present some of our results in Fig. ^ 





i(d) 

V 


n. = 0.0001 

. n. = 0.0005 

n =0.001 

n - 005 


^--"'^^tA-- 



0.005 0.01 0.015 20 40 60 80 100 

\x I A T(K) 

FIG. 4: (Color online) (a) In-plane DC conductivity in the 
bilayer as a function of the chemical potential (in units of 
A), (b) In-plane conductivity in the bilayer as a function 
of frequency at T = 300K and /x = (in units of A), (c) 
Perpendicular DC conductivity in the multilayer as a function 
of the chemical potential (in units of A), (d) Semi- log plot 
of the transport anisotropy {(t\\/(t±) in the multilayer as a 
function of the T (in K) at /x = 0. 

The DC conductivities are given by (Jdc = ctq ^(/^^ /^) 
where /i is the chemical potential. In the bilayer: 

Sdcii = 2 / 6^(^|fc2){lm[GSA(M,k)]lm[GgB(M,k)] 

+ Im[GAiB2(M,k)]lm[GA2Bi(/i,k)]} . (10) 

Upon taking the limit /i ^ 0, using the asymptotic ex- 
pressions for the self-energies in Eq. ©, the off-diagonal 
propagators drop out and the in-plane DC conductivity 
per plane acquires a universal minimum value given by: 



'"11," 



(3A)(eV/^), 



(11) 



independent of the impurity concentration. Hence, as in 
the case of the single layer [lO| , we find a universal con- 
ductivity minimum albeit with a different value (in the 
single layer one finds (Jmin = (4/7r)(e^//i)). This result 
shows that in these systems the in-plane conductivity is 
always of order e'^/h per plane and disorder independent. 
The frequency-dependence of the conductivity in the 
bilayer shows some structure. For the cleaner systems 
a Drude-like peak appears at cj = due to thermally 
excited carriers. The second peak dit co = t± is due to in- 
terband transitions. The perpendicular transport in the 



multilayer has the amazing property that close enough to 
half- filling the transport is enhanced by disorder, as can 
be seen in Fig. El (c). Also, as one can clearly see in Fig.EJ 
(d), the anisotropy ratio becomes very large as the Dirac 
point is approached, and the cleaner the system the larger 
the anisotropy. Using a similar asymptotic expression as 
in Eq. one expects that the anisotropy diverges as 
T~^/^ exactly at the Dirac point. For small, but finite 
values of the chemical potential, the anisotropy is still 
enormous but saturates at a finite value at T = 0. No- 
tice, however, that at high temperatures electron-phonon 
scattering (not discussed here) becomes important and 
can substantially modify the transport. 

Conclusions. In conclusion, we have presented results 
for the electronic properties of disordered graphene mul- 
tilayers showing that the behavior cannot be described 
by Landau's Fermi liquid theory of metals. The uncon- 
ventional behavior includes divergent self-energies near 
the Dirac point, the vanishing density of states on the A 
sublattice, the non-intuitive feature that disorder can in- 
crease the out-of-plane transport and the high anisotropy 
of the system near half-filling. These properties show 
that graphene multilayers are a new class of materials 
with an unusual metallic state |l6^. 
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